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Abstract 

In this paper, we study the problem of interpolating a continuous 
function at (n + 1) equally-spaced points in the interval [0,1], using 
shifts of a kernel on the (l/n)-spaced infinite grid. The archetypal 
example here is approximation using shifts of a Gaussian kernel. We 
present new results concerning interpolation of functions of exponen¬ 
tial type, in particular, polynomials on the integer grid as a step en 
route to solve the general interpolation problem. For the Gaussian 
kernel we introduce a new class of polynomials, closely related to the 
probabilistic Hermite polynomials and show that evaluations of the 
polynomials at the integer points provide the coefficients of the inter- 
polants. Finally we give a closed formula for the Gaussian interpolant 
of a continuous function on a uniform grid in the unit interval (assum¬ 
ing knowledge of the discrete moments of the Gaussian). 


1 Introduction 

In the mathematical literature pertaining to radial basis functions, there 
have been mainly two approaches to constructing interpolants with Gaussian 

*This work was supported by EPSRC Grant EP/H020071/1, the University of Leicester 
via study leave for the second author and Missouri State University through a generous 
travel grant for the second author. 
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kernels. The first involves interpolating a function / G C(M) on the h-spaced 
grid KL by an interpolant of the form 

y^^a^'ip{z - x/h), 

z&Z 

where 

exp( 


Analysis of this so-called cardinal approximation has been done in a series 
of papers of Baxter, Riemenschnieder and Sivakumar [H [HI EHl EQ] . 

The second concentrates on interpolating a continuous function on a hnite 
subset y of a compact interval (e.g [0,1]). Under this circumstance, the 
interpolant one seeks is of the form 

^ay'il){y - x). 

y&Y 

There are multidimensional set-ups for both approaches. Modern mathe¬ 
matical literature abounds in developing error estimates for approximation 
schemes in this context. We refer readers to [mOSlilHlE] and the references 
therein. 

Approximation methods involving sparse-grid algorithms have been re¬ 
cently proven effective and efficient; see 110]. Some sophisticated multi-level 
sparse grid kernel interpolation schemes have been constructed by authors 
of [H |9] . We are currently motivated to develop sparse-grid algorithms for 
high-dimensional approximation with the Gaussian kernel and derive error 
estimates for C^— functions with polynomial growth. However, there are 
several obstacles en route to achieving these goals. The main purpose of the 
current paper is to clear a few obstacles out of the way. First and foremost, 
we face the problem of interpolating a function at the {n + 1) equally-spaced 
points ih, i = 0,1, • • • , n with h = 1/n, where n G N. The approach we take 
here differs from those discussed in the above references. We hrst interpolate 
the given (n -|- 1) data by a degree n polynomial, and then interpolate the 
polynomial by a radial basis function interpolant on KL. 

Functions of exponential typ^ are often utilized as a half-way house in 
deriving error estimates for Sobolev space functions; see [13 ESI- As such, 

'^Some variants of functions of exponential type are also referred to as “band limited 
functions” in the literature. 
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it is worthwhile to study the effect of the interpolation scheme when the 
target functions are of exponential types, and in particular, polynomials, 
which we anticipate to play a significant role in our future effort to obtain 
more nuanced error estimates for functions with polynomial growth. 
Interestingly enough, we observe in the analytic number theory literature that 
interpolation goes the opposite way in the sense that functions of exponential 
type are employed to approximate the Gaussian and other useful radial basis 
functions; see EE]. We hope that interactions of the two seemingly inverse 
research tracks will create synergistic results. 

The layout of the paper is as follows. In Section |2l we will consider a 
general kernel V’ and study the operator T induced by the Toeplitz matrix 
V’(j ~ k), j,k & Z. The action of T on an / G C(M) takes the form: 

- k)f{j), /ceZ. 

j&i. 

We show that the operator is one-to-one on the space of polynomials, which 
are the subject of interest in the main body of the current paper. In the 
Appendix we show furthermore that the operator is one-to-one on a larger 
class of functions of exponential type. Moreover, we demonstrate that the 
coefficients of an interpolant can be written in closed form. In Section [21 we 
investigate the special case in which the Gaussian kernel is employed and the 
target functions are polynomials. We introduce new classes of polynomials 
resembling the classical probabilistic Hermite polynomials, and derive closed 
formulas for the coefficients of a Gaussian interpolant in terms of these poly¬ 
nomials. In Section SI we show how to use interpolate general functions on 
a uniform grid in the unit interval via interpolation by polynomials. This is 
not a convergent approximation scheme. However, it can be used as part of 
a residual correction scheme as is described in [9] . 

2 Interpolation with general kernels 

Let us assume we have an infinitely differentiable, positive function -0 such 
that all the “discrete moments” 

Mk = = ( 1 ) 

jez 
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are finite. Let Pk{.x) = for 0 < /c < iV. We seek coefficients akj, j G Z, 
such that 

I\Pk]{x) := ^ak,jij{x - j), (2) 

jei 

interpolates pk at all the integers. In other words 

/[p^](£)=£^, 0<A;<iV, ieZ. 

We will show that if coefficients (as polynomials in j of degree k) exist, 
then they are unique. We then construct such coefficients (of polynomial 
form). Thus I[pk] is the unique interpolant with coefficients of polynomial 
form. Furthermore, these coefficients are constructible in a recursive fashion. 

In the Appendix at the end of the paper we extend results of existence 
and uniqueness to include functions of exponential type. 

Lemma 1. Let be a polynomial of degree k for some Q < k < N. Suppose 
that 

f{i):=J2rkijmj-i)=0, ieZ. 

j&Z 

Then = 0. 

Proof: We proceed by induction. If /c = 0 then rk = c for some constant 
c. Then, 

-i) = 0 , 

jez 

and since fj is positive, c = 0. 

Assume that the result holds true for all polynomials of degree < k. Let 
Af{x) = f {x + 1) — f {x), X E M., he the forward difference operator. Then 

A/(£) = f{i+l)-f{i) 

= rk{j)ij{j - i) 

jez j&z 

= rkij + rk{j)f){j - i) 

j&Z j&Z 

= 5^(^fc(i + 1) - rk{j))f^{j - i) 
jez 

= '^Ark{j)f^{j - i) = 0. 

jez 
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Note that Ar^ is a polynomial of degree k — 1, which satishes Arfc(j) = 0, 
j G Z. By the induction hypothesis, we have that Ar^ = 0. It follows that 
is a constant. We use the induction hypothesis once more to conclude that 
Tfc is identically zero. □ 

The result of Lemma [T] can be equivalently stated as follows. There is no 
nontrivial polynomial p, deg(p) < N, such that 

J2p{mj-i) = 0, 

jez 

In the sequel, we will use the above fact without further declaration. 
Lemma 2. For /c = 0,1, • • • , N, 




i =0 


£ez. 


Proof: Substituting the expression for from ([T]) we have 


i=0 




ik-ir = 


i=0 


E , Ej‘-V(j) 


- jez 
k 




: ]3 


k—ini 


i=0 

j&Z 

jez 

by simply renumbering the sum. □ 

Let Bk,i = Mk-i{’lZl), k = 1- ■ ■ ,N + 1, i = 1, - ■ ■ ,k. Since 


k 

V(i-f) = E+ 

jez i=i 

we can use backward elimination to see that 

jez I i=i J 


( 3 ) 

( 4 ) 
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for some numbers j, fc = 1, • • • , A^+1, i = 1, • • • , fc. Define the polynomial 
ak{j) = Ak+i^i+if~^. Then the coefficients dehned in ([2]) satisfy 

An immediate consequence of Lemma [T] and (jl]) is 

Theorem 1. The coefficients Okj of the interpolant I\pk] given in ^ are 
unique. Moreover, they are polynomials (in j) of degree k — 1. 

Let 


Qkiffi = y^^Mk-m (5) 

m=0 '' 

In the next theorem we prove a generating function style relationship for the 
polynomials qk and the interpolating functions I[pk]- As a corollary of this 
theorem we give a recursive formula for the errors between I\pk] and pk- 

Theorem 2. For k = 0,1, ■ ■ ■ ,N, 

jez i=o ^ ' 

Proof: If we substitute dl]) into Lemma [2] we have 

jeZ i=0 ^ j€Z 

and by the linear independence of — £), £ E Z, 

f = ^(^.)Mk-Mj), jez. (7) 

i=0 
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We now expand the left hand side of ([S]) to get 


-xfipij - 


^'1 = EE 

j£Z i=0 


" ] {-xYf V(j - x) 


' k—i 


-^Y E 


^m=0 


k — i 


m 


EE 

jei, i=o 

E(:)(-)‘ECri^^‘--” 

1=0 ^ Y m=0 ^ 

k \ k—i 

0* X] ( _ * ]Mk-i-mI[Pm]{x), 


i—(j ) I 


- x) 


E 

i=0 


-X 


m=0 


m 


In deriving the above equations, we have used ([2]) and ([7]). Making use of 
the formula 


k\ fk — i 


m 


k\ fk — m 


m 


we reorder the hnal sum above as follows. 

k k—m 

- xf'ipij - x) = 


jez 


m=0 

k 


i=0 


k\ fk — i 


m 


Mk-m-i{-xy 


m=0 

k 


k—m 




i=0 


k — m 


Mk-m-i{-xY 


^ (^\l[Pm\{,x)qk-m{-x). 

rv)—n E / 


m=0 

In deriving the last equation, we have used (jS]). □ 

We dehne the error of interpolation Ek{x) ;= I\pk\{x) —pk{x), and 

Xk{x) := - x)’‘fj{j -x)-Mk, A; = 0,1, • • • , iV. 

jez 

Next, we obtain a recursive formula which, upon an appropriate inversion, 
allows us to write the errors in terms of the functions Xk and qk- 

Corollary 1. For k = 0,1, ■ ■ ■ ,N, 

k 


Xkix) = '^y.jqi{-x)Ek-i{x). 
i=0 '^*7 
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Proof: Since I[pk-i]{i) = Pk-i{i), i = 0,1, ■ ■ ■ , k, i E Z, hj Theorem HI 
we have for A: = 0, 1, • • • ,N, 

- f) = Mt. 

i=o ^ ' jez 

Since a polynomial is uniquely determined by its values on the integers we 
have 


= Mk. 

i=0 

Subtracting this equation from ([HD we see that 

gi(-x)Ilpk-i](x) - Mk 

= ^ (^^^(li{-x)I[Pk-i]{x) - (^'^^qi{-x)pk-i{x) 

i=0 ^ ' i=0 ^ ' 

qi{-x){I[pk-i\{x) -pk-i{x)) 

i=0 '' 

Remark 1. Theorem gives a recursive formula for computing the inter- 
polant of any polynomial, as long as one knows the interpolants for polyno¬ 
mials of lower degree. Likewise, the result of the corollary expresses the error 
between a polynomial and its interpolant in the same fashion. If the Gaus¬ 
sian kernel is employed to do interpolation, then we will have more interesting 
information to offer. We will study this topic in the next section. 

3 The Gaussian Kernel 

In this section we study exclusively the case of fj^x) = Pertinent 

to the contents of this paper will be the probabilistic Hermite polynomials 



- ? (■) 



Hck, k = 0,1, • • •. These may be defined in a number of ways, but for us 
perhaps the most appropriate one is via Rodrigues formula: 



We have the following explicit representation of these polynomials (see e.g. 


my- 


lk/2\ 


h — Oi 





where [xj is the greatest integer less than or equal to x. The polynomial Hck 
has a close cousin that is often referred to as the probabilistic polynomial of 
negative variance: 


Lfc/2J 


h — 97 



z=0 


which has the same coefficients in absolute value, but the coefficients are all 
positive. 

The probabilistic polynomials of negative variance arise very naturally in 
this study as they are the result of the continuous convolution of the Gaussian 
with the polynomials of appropriate degree: 

Lemma 3. Let'ip{x) =exp(—x^/2). Then 




where 



k even, 
k odd. 


Proof: It is well-known (see e.g. ira) that 



k even, 
k odd. 
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To prove the first equation, we make a simple change of variable w = y — x\ 

noo noo 

/ y^4>{x — y)dy = / {w + xYip{w)dw 


i=0 

k 

E 

j=0 

k 

E 

j=0 

\k/2\ 




Ck-ix" 


CiX 


k—i 


Eh..- 


^k—2i 


1=0 

We remind readers that all the odd degree terms have zero coefficients. If 
we substitute the value for C 2 i we see that Nek is the probabilistic Hermite 
polynomial of negative variance: 

~ *K^-2i)! 2* ■ 

A fascinating relationship between Hck and Nck is the so-called umbral 
composition (see 0 )- 


v--n\ 

E 


i=0 


k\ Hek-2i{x) 
i\{k-2i)\ 2 * 


Lfc/2J 

E(-i)‘. 

1 

2 = 0 


k\ Nek-2i{x) _ ^k 

. — (2/ . 


i\{k — 2i)\ 


Using Lemma [3] and the second equation above we have 


X 


k 


[k/2\ 

E(-i)‘ 


i=0 


k\ 

2H\{k-2i)\ 




Hek{y)'ilj{y 


x)dy, 


y)dy 


( 8 ) 


so that we can recover the monomials by integrating against the probabilistic 
Hermite polynomials. Of course, this gives us an idea of what will happen 
in the discrete case. 
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To do this, we need an analogue of the probabilistic Hermite polynomial 
for the discrete case. We dehne 


Hek{x) 



0 , 1 ,-, 


where M 2 k is the discrete moment as dehned in ([T]). We also let 


Nck = qk, /c = 0,1, • • • , 


(9) 


where qk is as dehned in Q. In other words 


Nek{x) = 



( 10 ) 


Equation ([8]) suggests that a closed formula for the interpolant I[pk]{x) 
resembles 

^Heki 3 )^/J{j - x), xeM. 
jez 

In the next result, we will show that this is indeed the case (up to a con¬ 
stant very close to 1) for k = 0, 1,2,3. For fc > 5 we need to make some 
modihcations to Hk for a closed form. 


Lemma 4. For k = 0,1, ■ ■ ■, and £ G Z, 

\k/i\ ( 2i 


'^Hek{z)'ilj{j - 


\k/i\ f 2i 

E B-D’ 

i=0 l^m =0 


k \ f k — 2m 


2m} \Ai — 2m 


M2mMii-2m ^ j 


k—Ai 


Proof: By Theorem [21 Equations ([2]) and fITUD . we have, for £ G Z, 




Lfc/2J 


^Hek{z)'ip{j - tj = 

jeZ [ i=o 
Lfc/2J 




E (-!)• (o) ^2, 15^ - <)1 

i=0 ^ ^ IjGZ J 



2 —2m 
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Rearranging we obtain 


- i) = 

jez 



where we have nsed the fact that, if i is odd, then 


£k—2i 


nk—Ai 


3=0 

which is trne becanse 
k \ f k — 2j 


‘^jJ - 2j 




k\ 


{k-2j)\ 


{k-2j)\{2j)\{k-2t)\{2t-2j)\ 

k\ {k — 2i + 2j)\ 

{k-2i + 2 j)\{ 2 i- 23 )\ {k- 2 i)\{ 23 )\ 


k 

2(*-j) 


k-2{i- j) 
2j 


□ 


As we can see, for A; = 0,1, 2, 3, the above lemma gives an exact formnla. 
Interestingly, snppose that we replace Mk by fc = 0,1, • • •, then the cor¬ 
rection terms above are all zero. This is why we get the nmbral composition 
formnla for the probabilistic Hermite polynomials. For higher degrees we 
need to modify the polynomial in the snmmation for interpolation. To this 
end we introdnce the polynomials Qk{x), which we dehne by 


Qk = A: = 0,1,2,3, (11) 

and for A: = 4, 5, • • •, 

[fc/4j f 2i / k \ / k 2 \ 1 

^ ~ ^ \2jJ J 

Using Lemma m we immediately get the main resnlt of this section 


Theorem 3. For A; = 0,1, • • •, and tp the Gaussian, we have 


i[pk]{x) = - x). 

jez 
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4 Gaussian Interpolant on a /i-spaced points 
in an interval 

In this section, we give a recipe for computation of a new Gaussian kernel in¬ 
terpolant to a function defined at equally spaced points X = {0, , 1}. 

The construction of this interpolant utilizes the full (1/n)— spaced inhnite 
grid. As such, it is different from most of Gaussian kernel interpolants con¬ 
structed with conventional procedures. However, for all the practical com¬ 
putational purposes, only a small number of centres outside of the interval 
of interpolation are required. The rapid decay of the Gaussian kernel offsets 
the error incurred by dropping terms (shifts) of the interpolant far from the 
interpolation interval. 

In this case we seek an interpolant of the form 

4[/](a:) = -nx), x e [0,1], 

jez 

for / G C([0,1]). To do this we follow the following recipe: 

1. Interpolate f on X with a degree n polynomial 

n 

Pn[f]{x) = '^Ci{f)ti{x), 
i=0 

where tj, i = 0, • • ■ , n is a basis for the degree n polynomials. 

2. Interpolate the scaled polynomial Si/nPn ■ x t—)■ P„(- x), at the integers 

0,1, • • • , n, and evaluate the result at nx: 

^n[f]{x) = I[Si/nPn]{nx) 

n 

= y^^Ci{Si/nf)I[ti]{nx). 
i=0 

With the basis of monomials this becomes 

n 

4[/](a^) = -nx), 

where Qi, i = O,--- ,n are as given in flTT]) and flT^ . The coefficients in 
the expression must be interpreted as the appropriate ones for the monomial 
basis. 
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5 Conclusion 


As main results of this paper, we have shown that the interpolant to a polyno¬ 
mial using a suitable kernel has polynomial coefficients. More importantly, 
a kernel interpolant to a polynomial is constructible recursively, as is the 
way in which we express the error between the polynomial and its kernel 
interpolant. For the Gaussian kernel, we provide closed formulas for the co¬ 
efficients of the kernel interpolant to a polynomial. These are given in terms 
of a new class of polynomials that closely resemble the classical probabilistic 
Hermite polynomials. Via interpolating polynomials, we hnd a way to con¬ 
struct a kernel interpolant to a function dehned on an equally-spaced grid of 
a compact interval. In theory, this interpolant uses shifts of the kernel on a 
full inhnite grid. In numerical implementation, however, only a small number 
of shifts of the kernel centered outside of the interpolation interval is needed 
thanks to the rapid decay of the kernel. These have cleared the way for our 
future work in which we will investigate numerical aspects of this process. 
Our goals are to obtain stable and efficient algorithms for the computation 
of the interpolant, and to develop error estimates for C'*'—functions having 
polynomial growth. With the stationary interpolation scheme, the error will 
not go to zero as the grid spacing contracts, but the errors estimate will be 
useful for analysing the residual approximation algorithm that is detailed in 
0 . 
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Appendix: Functions of exponential type 

In this appendix we show that we can extend the results of Section [2] 

beyond the polynomial case to that of functions of exponential type. For 

/ G T(M), we use the following Fourier transform pair: 



(13) 


We assume that both Fourier transform and inverse Fourier transform have 
been properly extended to the Schwartz class of tempered distributions. Let 
W denote the collection of all functions -0 G C(M) that decays rapidly. That 
is, there exists a constant G > 0, such that for any Y G N, the following 
inequality holds true: 



( 14 ) 
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Each G W induces a periodic function ip on R: 

:= xeR. (15) 

2:GZ 

The period of the above function is 1. We will use [0,1] as the fundamental 
interval. We are particularly interested in the subset W* of W dehned by 

W* := {V’ e W : pj{x) 7 ^ 0, X e [0,1]}. (16) 

Lemma 5. For each pj G W*, there exists a sequence of complex numbers 
Oz, z G li, such that J^zez ^ h G and 

X G R. 


'ip{x) 


= y^g^e 

2GZ 


2'Kizx 


Proof: First off, Wiener’s lemma [Tj p. 228] asserts that there exists a 
sequence of complex numbers Oz, z & Z, such that l^zl < 


1 

fj^x) 


a;GR. 

jsGZ 


The rapid decay of the function ip ensures that the periodic function pj (see 
([15])) is inhnitely differentiable on R. Since 'iplx) ^ 0 for all x G R, this prop¬ 
erty of smoothness of the function ip passes on to the function , whose 

Fourier coefficients a^, z G Z, therefore enjoy the desired decay condition. 

□ 

We will refer to the inequality as displayed in flT^ Wiener’s condition. Let 
S and S' denote, respectively, the Schwartz classes of functions and tempered 
distributions. For each given 0 < a < cxo, let Eg- denote the class of analytic 
functions of exponential type a. We will focus on a subclass E* of E^ defined 
by: 

E; := {/ G E, ; aC > 0 and at g M such that |/(x)| < C{1 + jx]^), x G R}. 

Proposition 4. Let ip G W*. For each a > 0, and every / G E* there exists 
a (7 G E*, such that 


fU) = '^ 9 iz)ip{j - z), jGZ. 
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Proof: By the Paley-Wiener theorem [2^ p. 162], we may write 
/ = f, T e5',supp(T) C [-a/(27r),a/(27r)]|l 


Since ^|J G W*, we have 

1 


xeR, 

zGTa 


in which the Fonrier coefficients Oz, z G Z, decay rapidly thanks to Lemma 
[5l That is, for any /c G Z+, we have \^zZ^\ < oo. Thus, the following 

equation dehnes ^ as a Schwartz class distribution: 

(y,0) := • 0), 0 G 5. 

Here denotes the function x i—)■ We also have 

supp C [-a/(27r),a/(27r)]. 

We use the same Paley-Wiener Theorem mentioned above to conclude that 
g := is an element of E* . Thus, there exists a constant C > 0 and an 
N eN such that 

g{x) < C{1 + |a;|^), x ^R. 

Now consider the function 


■= - z)'ip{z). 


Fix each hxed M > 0 and every x G [— M, M], we have 

^g{x-z)ij{z) 

zGlt 

<(75^(1 + 11-21") (1 + |2|("+2I)-‘ 

zez 

zez 

^The factor (27r)“^ is the result of the particular format of the Fourier transform pair 
we use in the present paper; see Equation m- 
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Here Cat is a constant depending only on N. Thns the series converges nni- 
formly on every compact snbset of M. Therefore the fnnction f* is continnous 
on R and has at most polynomial growth. We calcnlate its (distribntional) 
inverse Fonrier transform: 

(rr = (gY y = ^ y = r. 

This shows that both / and f* are the Fonrier transform of the distribntion 
T, meaning that they are the same fnnction. In particnlar, we have 

gU - z)ijiz) = g{z)'iij{j - z) = /(j), i e z. 


This completes the proof. □ 

For the nniqneness of the coefficients, we have the following result. 

Proposition 5. Assume that 0 < e < vr. Let g G '0 G W*, and let f 

be defined by 

fix) ■= Y9i^)f’ix - z), 

Then, in order that f{j) = 0, j G Z, it is necessary and sufficient that 
giO = 0, CeC. 

Proof: Of course, only the necessity part needs a proof. Assume that 
fU) = 0, j G Z. Write 

fid) = Y9iz)f^id -z) = Y^id - z)^fiz), j G Z, 

zGli zGli 

and consider the function f* dehned by, 

f*i^) = Y^i^ - z)f^iz)- 

Z&7. 

We remind readers that f{j) = f*ij), j G Z. From the proof of Proposition 
m we observe that f* is continuous on R and has at most polynomial growth. 
The (distributional) Fourier transform of f* can be easily calculated to be 

/•«) - sK) = m-m- 

zez 
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Since g G E*_^, g is supported in [—1/2 + e, 1/2 — e]. Thus f* G E*_^. 
Resorting to Carlson’s theorem]^ we conclude that f*{C) = 0) C ^ C. The 
Fourier transform of f* is therefore also zero. Since /*(C) = and 

V’(C) 7^ 0) we have 'g{C) =0, C ^ That is, 'g is the zero distribution. Hence 
g is identically zero. □ 

Propositions m and [5] imply the following result. 

Corollary 2. Let 0 < e < tv be given, and let ip G W*. For each g G 
there exists a unique f G such that 

= gU), j e z. 

z£Z 

Suppose that / is radial (even), and that for some 5 > 0 we have 

\f{x)\ < A(1 + |7K)| < ^(1 + 

where H > 0 is a constant. Then the following Poisson summation formula 
holds true; see [2T1 p.252]. 

zGX zGIa 


Thus, Wiener’s condition is satished if both -i/ and ip have the decay rate 
shown in flTT)) . and Tp is positive. Specihcally, the Gaussian kernel satishes 
this condition. 

We remind readers that for any 0 < e < vr, and any polynomial p, we 
have p G E*_g. To interpolate a polynomial, we do not need any extra 
decay condition other than what has been imposed on functions from W*. 
Moreover, any -0 G iS satisfying Wiener’s condition suffices. 


•^Carlson’s theorem asserts that a function in that vanishes on all the positive 

integers is identically zero; see [18]. 
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